rm(list=ls(all=TRUE))  # new
library(foreign)
library(pscl)

load("rc.Rdata")
load("results.Rdata")
ideal.est = results[[1]]




### make a matrix with information on each legislator

trim <- function (x) gsub("^\\s+|\\s+$", "", x) # returns string w/o leading or trailing whitespace

make.members.matrix <- function(){
	members1 = colnames(ideal.est$x)
	members2 = strsplit(members1,split="[.]")
	for (i in 1:length(members2)){
		m = members2[[i]]
		mname = trim(paste(m[5:length(m)],collapse="."))
		members2[[i]] = c(m[1:4],mname)
	}
	t(matrix(unlist(members2),nrow=5))
}
members = make.members.matrix()
colnames(members)=c("type","cong","id","party","name")

### positions of members satisfying the criteria
get.pos <- function(type="",cong="",id="",party=""){
	cong1 = format(cong)
	id1 = format(id)
	pos1 = (members[,"type"] %in% type) | type==""
	pos2 = (members[,"cong"] %in% cong1) | cong==""
	pos3 = (members[,"id"] %in% id1) | id==""
	pos4 = (members[,"party"] %in% party) | party==""

	which(pos1&pos2&pos3&pos4)
}

### simulated ideal points of members satisfying the criteria
get.x <- function(type="",cong="",id="",party="",print.members=TRUE,save.att=FALSE){
	pos = get.pos(type,cong,id,party)
	att = members[pos,]
	if (print.members) {
		if (length(pos)==1) {
			print(cbind(pos,t(att)))
		} else { 
			print(cbind(pos,att))
		}
	}
	if (save.att){
		if (length(pos)==1){
			return(rbind(cong = att["cong"],type = att["type"],
			 party=att["party"],ideal = ideal.est$x[,pos,1]))

		} else {		
		return(rbind(cong = att[,"cong"],type = att[,"type"],
		 party=att[,"party"],ideal = ideal.est$x[,pos,1]))
		}
	} else{
		return(ideal.est$x[,pos,1])
	}
}


### search by name (partial match)
search.by.name <- function(mname){
	pos = grep(mname,members[,"name"])
	cbind(pos,members[pos,])
}


#### EXAMPLE ###################################################
#search.by.name("S[EU]T")  # search for SET or SUT
#search.by.name("S..TT")  # search for S??TT
#search.by.name("OSTETT")  # search for OSTETT
################################################################
#x1=get.x(id="29520")
#x105=get.x(id="29520",cong=105)
#x106=get.x(id="29520",cong=106)
#xHR103=get.x(cong=103,party="R",type="H")
#xP=get.x(type="P")
#x105_106 = get.x(id="29520", cong=c(105,106))
########################################################################


### posterior prob. that person "a" is more liberal than person "b" (A and B might 
###be the same person in different Congresses)

prob.more.lib <- function (a.id, a.cong, b.id, b.cong){
	a = get.x (id=a.id, cong=a.cong, print.members=FALSE)
	b = get.x (id=b.id, cong=b.cong, print.members=FALSE)
 	mean (a<b) 
}


##################################
# RUN
##################################

########################################
# Create matrices with members' ideal points for each congress
#
########################################
#dw.mat

dw = read.dta(file="SL01112D21_BSSE.DTA") #DW-NOMINATE for the Senate

# col1: congress num
# col2: id 
# col4: Congressional District Number (0 if Senate or President)
# col6: party (100=D, 200=R)
# col8: 1st dim coordinate

subset = dw [dw[,1]>100 ,]
dw.id.table = rownames(table(subset[,2]))

dw.mat = matrix(NA, nrow = length(dw.id.table), ncol = 14)
rownames(dw.mat) = dw.id.table
colnames(dw.mat) = c( 101:112, "party","type")

for (j in 1:length(dw.id.table)){
	for (i in 1:12){
		id.pos = subset[,2]%in%dw.id.table[j]
		pos.list = which(id.pos)
		
		for (k in 1:length(pos.list)){
	#browser()	
		cong = subset[pos.list[k],1]
		
		dw.mat [j,(cong-100)]= subset[pos.list[k],8]
		dw.mat [j,"type"]= subset[pos.list[k],4]
		dw.mat [j,"party"]= subset[pos.list[k],6]
		
		}
	}
}

## Standard errors!!
dw.mat.se = matrix(NA, nrow = length(dw.id.table), ncol = 14)
rownames(dw.mat.se) = dw.id.table
colnames(dw.mat.se) = c( 101:112, "party","type")

for (j in 1:length(dw.id.table)){
	for (i in 1:12){
		id.pos = subset[,2]%in%dw.id.table[j]
		pos.list = which(id.pos)
		
		for (k in 1:length(pos.list)){
	#browser()	
		cong = subset[pos.list[k],1]
		
		dw.mat.se [j,(cong-100)]= subset[pos.list[k],"bootse1"]
		dw.mat.se [j,"type"]= subset[pos.list[k],4]
		dw.mat.se [j,"party"]= subset[pos.list[k],6]
		
		}
	}
}



# change.mat

###################################################################
## change.mat is a matrix with posterior medians. 
## Each column from the third is for each congress, and each row is for each id. 
## For example, change.mat[1,3] shows the posterior median of 
#the first legislator in the 101st congress.
###################################################################
change.mat = matrix(NA, nrow = length(dw.id.table), ncol = 14)
rownames(change.mat)=dw.id.table
colnames(change.mat)=c("type","party",101:112)


for (i in 1:length(dw.id.table)){
	temp = get.x(id=dw.id.table[i],save.att=TRUE, print.members=FALSE)
	n.cong = dim(temp)[2]

	for(j in 1:n.cong){
		congr = as.numeric(temp["cong",j])-100
		change.mat[i,1]=temp["type",j]
		change.mat[i,2]=temp["party",j]
		change.mat[i,(congr+2)] = median(as.numeric(temp[c(-1,-2,-3),j]),na.rm = TRUE)	
	}
}	


#numbers only
dw.mat.num = dw.mat[,c(-13,-14)]
change.mat.num = change.mat[,c(-1,-2)]

### Party Polarization

h.dem.vec = (change.mat[,"type"]=="H" & change.mat[,"party"]=="D")
h.rep.vec = (change.mat[,"type"]=="H" & change.mat[,"party"]=="R")
s.dem.vec = (change.mat[,"type"]=="S" & change.mat[,"party"]=="D")
s.rep.vec = (change.mat[,"type"]=="S" & change.mat[,"party"]=="R")

# dw_nominate
h.dem.med.dw = rep(NA, times=12)
h.rep.med.dw = rep(NA, times=12)
s.dem.med.dw = rep(NA, times=12)
s.rep.med.dw = rep(NA, times=12)

for (i in 1:length(h.dem.med.dw)){

	h.dem.med.dw[i] = median(dw.mat[h.dem.vec,i],na.rm=TRUE)
	h.rep.med.dw[i] = median(dw.mat[h.rep.vec,i],na.rm=TRUE)
	s.dem.med.dw[i] = median(dw.mat[s.dem.vec,i],na.rm=TRUE)
	s.rep.med.dw[i] = median(dw.mat[s.rep.vec,i],na.rm=TRUE)
}

#ours
h.dem.med = rep(NA, times=12)
h.rep.med = rep(NA, times=12)
s.dem.med = rep(NA, times=12)
s.rep.med = rep(NA, times=12)

for (i in 1:length(h.dem.med)){

	h.dem.med[i] = median(as.numeric(change.mat.num[h.dem.vec,i]),na.rm=TRUE)
	h.rep.med[i] = median(as.numeric(change.mat.num[h.rep.vec,i]),na.rm=TRUE)
	s.dem.med[i] = median(as.numeric(change.mat.num[s.dem.vec,i]),na.rm=TRUE)
	s.rep.med[i] = median(as.numeric(change.mat.num[s.rep.vec,i]),na.rm=TRUE)
}



###############################
#
#  FIGURE 8
# Presidential Candidate
#
###############################

par(mfrow=c(1,2))
#### DW-NOMINATE
#(1) McCain
mccain.id = search.by.name("MCCAIN")[,"id"][1]
mccain = dw.mat[(rownames(dw.mat)==mccain.id),]
mccain.se = dw.mat.se[(rownames(dw.mat)==mccain.id),]


plot(x=102:111, s.rep.med.dw[2:11],lwd=1, type="l",col=1,ylim=c(-0.8,0.7),
pch=15,xlab="Congress",ylab="Ideal point estimate", main=expression("DW-NOMINATE"))
points(x=102:111, y=mccain[2:11],lwd=2,pch=c(19,19,19,19,1,19,19,19,1,19))
legend("topleft", legend=c("McCain","McCain running for president"), pch=c(19,1), bty="n")
legend("topright", legend=c("Sen Republican Median"),lty=1,lwd=1, box.col=0)

arrows(102:111, mccain[2:11]-1.96*mccain.se[2:11], 102:111, mccain[2:11]+1.96*mccain.se[2:11],length=0.05, angle=90, code=3)


#(2) KERRY
kerry.id = search.by.name("KERRY")[,"id"][1]
kerry=dw.mat[(rownames(dw.mat)==kerry.id),]
kerry.se = dw.mat.se[(rownames(dw.mat)==kerry.id),]
 
lines(x=102:111, s.dem.med.dw[2:11], lty=2,lwd=1)
points(x=102:111, y=kerry[2:11],lwd=2,pch=c(17,17,17,17,2,17,17,17,2,17))
legend("bottomleft", legend=c("Kerry running for president","Kerry"), pch=c(2,17), bty="n")
legend("bottomright", legend=c("Sen Democrat Median"), bty="n",lty=2,lwd=1)
arrows(102:111, kerry[2:11]-1.96*kerry.se[2:11], 102:111, kerry[2:11]+1.96*kerry.se[2:11], length=0.05, angle=90, code=3)

#### Anchors Away

# (1) McCain
McCain.info = search.by.name("MCCAIN")
McCain = get.x(id= rownames( table(McCain.info[,"id"]) ) )


#
mccain.sub =McCain[,2:11]
colnames(mccain.sub)=c(102:111)
boxplot(mccain.sub, col=c(0,0,0,0,1,0,0,0,1,0),notch=TRUE, xlab="Congress", ylim=c(-2.2,2),
ylab = "Ideal point estimate", main=expression(paste(plain('"Anchors Away"'),plain('ideal point estimate'))))
points(s.rep.med[2:11],type="l",col=1,lwd=1)
legend("topleft", legend=c("McCain","McCain running for president"), 
fill=c(0,1), box.col=0)
legend("topright", col=1, legend=c("Sen Republican Median"), box.col=0, lty=1, lwd=1)



# (2) Kerry
Kerry.info = search.by.name ("KERRY")
Kerry = get.x(id = rownames( table(Kerry.info[,"id"]) ) )

Kerry.sub =Kerry[,2:11]
colnames(Kerry.sub)=c(102:111)	
boxplot(Kerry.sub, notch=TRUE,col=gray(c(0.5,0.5,0.5,0.5,0.5,0.5,0.9,0.5,0.5,0.5)),add=TRUE)
points(s.dem.med[2:11],type="l",lty=2,lwd=1)


legend("bottomleft", legend=c("Kerry","Kerry running for president"), 
fill=gray(c(0.5,0.9)), box.col=0)

legend(x=0,y=-1.6, legend=c("Sen Democrat Median"), box.col=0, lty=2, lwd=1)



